#include <vector> 
#include <map>
#include <list>
#include <stdio.h>
#include <string.h>
#include <stdio.h>

#include <cstdlib>
#include <string>
#include <cmath>
#include <iostream>
#include <fstream>
#include <sstream>
#include "ulib.h"
#include "graph.h"

#define WARP 16
#define MAXLOAD 8

using namespace std;


void init();
int bc_vertex(int *h_ptrs, int* h_js, int n_count, int e_count, int nb, float *h_bc);
int bc_edge(int* h_v, int *h_e, int n_count, int e_count, int nb, float *h_bc);

template <typename VtxType, typename WeightType>
void ReadGraphFromFile(FILE *fpin, VtxType *numofvertex, VtxType **pxadj, VtxType **padjncy, WeightType **padjncyw, WeightType **ppvw);
template <typename VtxType, typename WeightType>
void ReadGraphFromMMFile(FILE *matfp, VtxType *numofvertex, VtxType **pxadj, VtxType **padjncy, WeightType **padjncyw, WeightType **ppvw);
static int really_read(std::istream& is, char* buf, size_t global_size);

template <typename VtxType, typename EdgeType, typename WeightType>
void ReadBinary(char *filename, VtxType *numofvertex_r, VtxType *numofvertex_c, EdgeType **pxadj, VtxType **padj, WeightType **padjw, WeightType **ppvw);
template <typename VtxType, typename WeightType>
void ReadGraph(char *filename, VtxType *numofvertex, VtxType **pxadj, VtxType **padjncy, int** ptadj, WeightType **padjncyw, WeightType **ppvw);


void bc_cpu (int* xadj, int* adj, int nVtx, int nz, int nb, float* bc) {

	for (int i = 0; i < nVtx; i++)
		bc[i] = 0.;

	int* bfsorder = new int[nVtx];
	int* Pred = new int[xadj[nVtx]];
	int* endpred = new int[nVtx];
	int* level = new int[nVtx];
	int* sigma = new int[nVtx];
	float* delta = new float[nVtx];

	for (int source = 0; source < min (nb, nVtx); source++) {
		int endofbfsorder = 1;
		bfsorder[0] = source;

		for (int i = 0; i < nVtx; i++)
			endpred[i] = xadj[i];

		for (int i = 0; i < nVtx; i++)
			level[i] = -2;
		level[source] = 0;

		for (int i = 0; i < nVtx; i++)
			sigma[i] = 0;
		sigma[source] = 1;

		//step 1: build shortest path graph
		int cur = 0;
		while (cur != endofbfsorder) {
			int v = bfsorder[cur];
			for (int j = xadj[v]; j < xadj[v+1]; j++) {
				int w = adj[j];
				if (level[w] < 0) {
					level[w] = level[v]+1;
					bfsorder[endofbfsorder++] = w;
				}
				if (level[w] == level[v]+1) {
					sigma[w] += sigma[v];
				}
				else if (level[w] == level[v] - 1) {
					Pred[endpred[v]++] = w;
				}
			}
			cur++;
		}

		for (int i = 0; i < nVtx; i++) {
			delta[i] = 0.;
		}

		//step 2: compute betweenness
		for (int i = endofbfsorder - 1; i > 0; i--) {
			int w = bfsorder[i];
			for (int j = xadj[w]; j < endpred[w]; j++) {
				int v = Pred[j];
				delta[v] += (sigma[v] * (1 + delta[w])) / sigma[w];
			}
			bc[w] += delta[w];
		}
	}

	delete[] bfsorder;
	delete[] Pred;
	delete[] level;
	delete[] sigma;
	delete[] delta;
	delete[] endpred;
}


int intcmp(const void *v1, const void *v2) {
	return (*(int *)v1 - *(int *)v2);
}


int intcmprev(const void *v1, const void *v2) {
	return (*(int *)v2 - *(int *)v1);
}


void compute_deg_dist (int* xadj, int n) {
	int* degs = (int*) calloc(n, sizeof(int));
	for (int i = 0; i < n; i++) {
		degs[xadj[i+1] - xadj[i]]++;
	}

	for (int i = 0; i < n; i++) {
		if (degs[i] > 0)
			cout<<i<<" "<<degs[i]<<endl;
	}

	free(degs);
	exit(1);
}


int main(int argc, char** argv) {

	init();
	char c;
	char* infilename, outfilename;
	int threads_per_block = -1, paropt = -1, nb = 1, xpar = 0, times = 1, sortopt = 0;

	if (argc != 7) {
		std::cout<<"usage: "<<argv[0]<<" <will_order:0|1> <deg1:0|1> <kernel> <nbsource> <xpar>"<<std::endl;
		return -1;
	}

	int will_order = atoi(argv[2]);
	int deg1 = atoi(argv[3]);
	int btype = atoi(argv[4]);

	nb = atoi(argv[5]);
	xpar = atoi(argv[6]); // max-deg # of a virtual vertex

	int n, i, j, nVtx;
	int *xadj, *adj, *tadj, *mark, *queue, *new_xadj, *new_adj;

	int count, maxcompid, u, v, p;
	int nocomp;
	int qptr, qeptr, largestSize, compsize, lcompid;

	ReadGraph<int, int>(argv[1], &n, &xadj, &adj, &tadj, NULL, NULL);
	nVtx = n;
	int* compid = (int*) malloc(sizeof(int) * n);
	for(i = 0; i < n; i++)
		compid[i] = -1;
	int* que = (int*) malloc(sizeof(int) * n);

	printf("there are %d vertices %d edges\n", n, xadj[n]);

	nocomp = qptr = qeptr = largestSize = 0;
	for (int i = 0; i < n; i++) {

		if(compid[i] == -1) {
			compsize = 1;
			compid[i] = nocomp;
			que[qptr++] = i;

			while(qeptr < qptr) {
				u = que[qeptr++];
				for(p = xadj[u]; p < xadj[u+1]; p++) {
					v = adj[p];
					if(compid[v] == -1) {
						compid[v] = nocomp;
						que[qptr++] = v;
						compsize++;
					}
				}
			}
			if(largestSize < compsize) {
				lcompid = nocomp;
				largestSize = compsize;
			}
			nocomp++;
		}
	}

	int nz = xadj[n];
	int ecount = 0;
	int vcount = 0;

	for(i = 0; i < n; i++) {
		if(compid[i] == lcompid) {
			que[i] = vcount++;
			for(p = xadj[i]; p < xadj[i+1]; p++) {
				if(compid[adj[p]] == lcompid)
					ecount++;
			}
		}
	}

	int* lxadj = (int*) malloc(sizeof(int) * (vcount+1));
	int* ladj = (int*) malloc(sizeof(int) * (ecount));
	int* ltadj = (int*) malloc(sizeof(int) * (ecount));
	vcount = 0; ecount = 0;
	lxadj[0] = 0;
	for(i = 0; i < n; i++) {
		if(compid[i] == lcompid)  {
			vcount++;

			for(p = xadj[i]; p < xadj[i+1]; p++) {
				if(compid[adj[p]] == lcompid) {
					ladj[ecount++] = que[adj[p]];
				}
			}
			lxadj[vcount] = ecount;
		}
	}
	printf("largest component graph obtained with %d vertices %d edges -- %d\n", vcount, ecount, lxadj[vcount]);


	n = vcount;
	nz = ecount;
	free(xadj); xadj = lxadj;
	free(adj); adj = ladj;
	free(tadj); tadj = ltadj;

	int* degs = (int*)malloc(sizeof(int) * n);
	int* myedges = (int*)malloc(sizeof(int) * nz);
	memcpy(degs, xadj, sizeof(int) * n);

	int ptr;
	for(i = 0; i < n; i++) {
		for(ptr = xadj[i]; ptr < xadj[i+1]; ptr++) {
			j = adj[ptr];
			myedges[degs[j]++] = i;
		}
	}

	for(i = 0; i < n; i++) {
		if(xadj[i+1] != degs[i]) {
			printf("something is wrong i %d xadj[i+1] %d degs[i] %d\n", i, xadj[i+1], degs[i]);
			exit(1);
		}
	}

	memcpy(adj, myedges, sizeof(int) * xadj[n]);
	for(i = 0; i < n; i++) {
		for(ptr = xadj[i]+1; ptr < xadj[i+1]; ptr++) {
			if(adj[ptr] <= adj[ptr-1]) {
				printf("is not sorted\n");
				exit(1);
			}
		}
	}

	memcpy(degs, xadj, sizeof(int) * n);
	for(i = 0; i < n; i++) {
		for(ptr = xadj[i]; ptr < xadj[i+1]; ptr++) {
			j = adj[ptr];
			if(i < j) {
				tadj[ptr] = degs[j];
				tadj[degs[j]++] = ptr;
			}
		}
	}

	free(degs);
	free(myedges);

	for(i = 0; i < n; i++) {
		for(ptr = xadj[i]; ptr < xadj[i+1]; ptr++) {
			j = adj[ptr];
			if((adj[tadj[ptr]] != i) || (tadj[ptr] < xadj[j]) || (tadj[ptr] >= xadj[j+1])) {
				printf("error i %d j %d ptr %d\n", i, j, ptr);
				printf("error  xadj[j] %d  xadj[j+1] %d\n",  xadj[j], xadj[j+1]);
				printf("error tadj[ptr] %d\n", tadj[ptr]);
				printf("error adj[tadj[ptr]] %d\n", adj[tadj[ptr]]);
				exit(1);
			}
		}
	}


	int* map_for_order = (int *) malloc(n * sizeof(int));
	int* reverse_map_for_order = (int *) malloc(n * sizeof(int));
	int* weight = (int *) malloc(sizeof(int) * n);
	float* bc  = (float *) malloc(sizeof(float) * n);

	for(i = 0; i < n; i++) {
		weight[i] = 1;
		map_for_order[i] = -1;
		reverse_map_for_order[i] = -1;
	}


	FILE* ofp;
	ofp = fopen("bc_out.txt", "w");

	nb = min(n, nb);

	printf("will be executed on %d vertices %d %d edges\n", n, xadj[n], nz); fflush(0);

	// kernels..
		if (btype == 0) {
			printf("BC: CPU based\n");
			bc_cpu (xadj, adj, n, nz, nb, bc);
		}
		else if (btype == 1) {
			printf("BC: vertex based parallelism\n");
			bc_vertex (xadj, adj, n, nz, nb, bc);
		}
		else if(btype == 2) {
			printf("BC: edge based parallelism\n");

			int* is = (int*) malloc(sizeof(int) * nz);
			for(i = 0; i < n; i++) {
				for(ptr = xadj[i]; ptr < xadj[i+1]; ptr++) {
					is[ptr] = i;
				}
			}

			bc_edge (is, adj, n, nz, nb, bc);

			free (is);
		}
		
		for (int i = 0; i < n; i++) {
				fprintf(ofp, "bc[%d]: %lf\n", i, bc[i]);
		}

	fclose(ofp);
	FILE* lfp;
	lfp = fopen("bc_out.txt", "r");
	char a,b,d,e,f,g, h;
	double val;
	int id;
	double* result_bc = new double[nVtx];
	while (fscanf(lfp, "%c%c%c%d%c%c%c%lf%c", &a, &b, &f, &id, &d, &e, &g, &val, &h) != EOF) {
	  result_bc[id] = val;
	}
	fclose(lfp);
	ofp = fopen("bc_out.txt", "w");
	for(int i = 0; i < nVtx; i++)
		fprintf(ofp, "bc[%d]: %lf\n", i, result_bc[i]);
	fclose(ofp);
	delete[] result_bc;

	free (bc);


	return 0;
}

